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ABSTRACT 


High resolution windspeed profile measurements are 
needed to provide reliable detection of hazardous low- 
altitude windshear with an airborne pulse Doppler radar, the 
system phase noise in a Doppler weather radar may degrade the 
spectrum moment estimation quality and the clutter 
cancellation capability which are important in windshear 
detection. Also the bias due to weather return Doppler 
spectrum skewness may cause large errors in pulse pair 
spectral parameter estimates. These effects are analyzed for 
the improvement of an airborne Doppler weather radar signal 
processing design. This dissertation also presents a method 
for the direct measurement of windspeed gradient using low 
pulse repitition frequency (PRF) radar. This spatial gradient 
is essential in obtaining the windshear hazard index. As an 
alternative the modified Prony method is suggested as a 
spectrum mode estimator for both the clutter and weather 
signal. Estimation of Doppler spectrum modes may provide the 
desired windshear hazard information without the need of any 
preliminary processing requirement such as clutter filtering. 
The results obtained by processing a NASA simulation model 
output support consideration of mode identification as one 
component of a windshear detection algorithm. 
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CHAPTER I 


INTRODUCTION 
Doppler Weather Radar 

Radar is considered to have a great potential as a 
remote sensing device, but ordinary radar can not see what 
the winds are doing. Only radars using the Doppler frequency 
shift can show the wind speed and its direction. These radars 
are called Doppler radars. Basically there are two different 
types of Doppler radars, i.e., a continuous wave (CW) Doppler 
radar and a pulse Doppler radar. The pulse Doppler radar has 
an advantage over the CW radar in that the detection 
performance is not limited by transmitter leakage or by 
signals reflected from nearby clutter or from the radome . The 
pulse Doppler radar avoids this difficulty since its receiver 
is turned off during transmission whereas the CW radar 
receiver is always on. Previous application of pulse Doppler 
radar techniques in mapping severe storm reflectivity and 
velocity structure has been very successful [1],[2]. Since 
Doppler radars have the capability of seeing dynamic 
structure of air in low reflectivity conditions [3], these 
radars are considered useful in the investigation of dry 
weather situations. 

A pulse Doppler weather radar transmits a train of 
pulses and receives reflected signals as in phase (I) and 
quadrature phase (Q) components to differentiate receding or 
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approaching targets. Figure 1 shows a simplified Doppler 
radar diagram and Figure 2 illustrates the transmitted 
waveform and the demodulated return signal. Moving targets in 
the turbulent air illuminated by the radar will shift the 
transmitted frequency which is well known as the Doppler 
effect. This Doppler frequency shift which represents the 
time rate of phase change in the demodulated signal is 
determined as 



where V r is the radial velocity of target particles and X. is 
the wavelength of transmitted signal. As seen in Figure 2, 
the change in signal phase typically is extremely small 
within a pulse duration (e.g., T= 1 (lsec, and weather target 

velocities on the order of tens of m/sec) . Hence, target 
phase shifts are measured over an interpulse period T s which 

means that the pulse Doppler radar behaves as a phase 
sampling device. 

From the I,Q data, important information (i.e., mean 
velocity, spectrum width, reflectivity, etc.) can be 
extracted using the appropriate radar signal processing 
techniques. In coherent Doppler radar systems, the accuracy 
of such information depends on the system phase stability. By 
coherent it is meant that the phase of transmitted signal is 
preserved in the reference signal. As seen from Figure 1, 



PULSE 



Figure 1. Simplified Block Diagram of a Doppler Radar 
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usually two highly stable oscillators, i.e, stable local 
oscillator (STALO) and coherent oscillator (COHO) , are needed 
in the heterodyne detection system to retain the phase 
coherence. Both oscillators are used to provide a transmitted 
and reference signal, but the stability problem of the COHO 
is not so serious compared with that of the STALO since the 
frequency range of the COHO is several orders of magnitude 
less than that of the STALO. Therefore it is considered that 
the STALO plays a major role in the pulse-to-pulse system 
phase stability problem. Of course, there may be many other 
factors increasing the system phase noise such as the 
mismatch problem in obtaining I, Q data [4]. Phase 
instability factors will degrade Doppler frequency resolution 
and the spectrum moment estimation quality. This phase noise 
effect may also severely limit the capability of clutter 
rejection filtering by spreading the clutter signal over the 
Doppler signal bandwidth as seen in Figure 3. 

Although Doppler radar operation was successful in some 
other experimental cases [5], [6], future application of 
Doppler weather radar where higher spatial resolution is 
needed may require higher Doppler resolution and more 
accurate estimation of spectral moments [50], [51]. A 
fundamental problem is quantification of the effect of radar 
system phase noise on the performance of Doppler weather 
radar. However, it is really difficult to investigate all 
phase noise contributing factors separately in a radar 
system. Thus, a general systematic approach was developed in 



clutter return 



a return spectrum without phase noise 



a return spectrum with phase noise 


Figure 3. Comparison of Two Doppler Return 

Spectra Showing Phase Noise Effect 
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this dissertation to analyze the whole system phase noise 
effect . 

Two popular methods used in radar signal processing for 
spectral moment estimation are the pulse pair method and the 
Discrete Fourier Transform (DFT) method. Modern spectral 
estimation techniques are considered to be much better than 
classical methods, but they may not be suitable for real time 
processing. The need to provide computed output as data are 
obtained from the radar may be essential in some applications 
such as airborne windshear detection when the look ahead time 
is limited by the ranging capability of the radar and may be 
on the order of tens of seconds at most . There are several 
methods which may be used in the estimation of weather 
spectral moments [7] . Among them, the most economical is 
generally considered to be the pulse pair estimator, often 
called a covariance estimator [8] . Especially in a weather 
radar system, this pulse pair algorithm is more widely 
accepted since it is simple to implement and fast enough to 
process huge amounts of data for real time mapping of the 
weather situation in an interested area. It is also shown in 
[52] that the performance of the pulse pair estimator is even 
better than that of the DFT estimator at low signal-to-noise 
(SNR) ratios and narrower widths. However, the DFT estimator 
has some advantages including the absence of bias due to 
nonsymmetric spectra and the feasibility of eliminating 
anomalous spectral powers , Also DFT based methods experience 
finite parameter estimation errors at very large spectral 
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widths which is not so with pulse pair processing where 
errors will increase exponentially as spectral widths 
increase. The other weakness in pulse pair estimation is that 
the algorithm may yield meaningless results in the case of a 
multimodal return spectrum while, for example, DFT processing 
will show all the modes in a return spectrum. In this 
dissertation these two popular estimators are reevaluated 
considering the system phase noise effect. Also, bias in the 
pulse pair estimation, introduced by a skewed spectrum, is 
analyzed. Alternatives for spectral parameter estimation 
which may be more robust in the presence of multi-modal 
spectra are also considered. 

Windshear 

One of the important potential applications of Doppler 
weather radar is in a windshear detection system. When the 
wind abruptly shifts its speed or direction, it can mean 
deadly difficulty for an airliner particularly at low 
altitudes such as on approach or take off. This dangerous 
windshear is frequently caused by microbursts. The term 
"microburst" was first used by Fujita [10] to describe a 
relatively small column of downward-rushing air when he was 
investigating the 1975 crash of an Eastern Airlines Boeing 
727 at New York's Kennedy International Airport. Since then, 
an estimated 26 major aircraft accidents between 1964 and 
1985 have been attributed to the microburst. Microbursts are 
sudden downdrafts of highly turbulent air that may cause very 
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hazardous windshear conditions. In their presence a plane can 
first encounter a sharp head wind, then an intense down 
draft, and finally a strong tail wind, all in a matter of 
seconds (see Figure 4) . These resulting wind currents appear 
as if they are designed to cause airline crashes. To make 
matters worse, a recent study by Fujita indicates that 
microbursts are more common than previously thought and can 
be created by relatively small harmless-looking rain clouds, 
not just large thunderstorms. Therefore, it is really 
difficult to avoid these dangerous situations without any 
early warning system designed to detect windshear. However, 
conventional weather radars may not be appropriate for use in 
this kind of situation. 

Since microbursts can occur within a very small 
geographical scale and the reflectivity of dry microbursts 
may be very weak, the weather radar for microburst detection, 
should have high sensitivity and high resolution of both 
range and Doppler frequency. Currently NASA is developing an 
airborne Doppler weather radar for the detection of 
microbursts, the number one killer of U.S. airline 
passengers. For the same purpose, the terminal Doppler 
weather radar (TDWR) is also now under development by 
Raytheon Co. for the Federal Aviation Administration (FAA) . 
The second major function of TDWR is to improve air traffic 
management through forecasts of windshifts, precipitation and 
other weather related hazards. However, they need further 
improvement and verification to be used as a windshear 
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detection system since undetected microbursts may cause 
tragic accidents. The other weather radar called N EXRAD (next 
generation weather radar), which is just beginning to be 
deployed nationwide, may have the microburst— induced 
windshear detection capability, but its main purpose is to 
replace the non-Doppler meteorological radars of the National 
Weather Service. NEXRAD will operate in the 5-6 GHz frequency 
range and use a low pulse repetition frequency (PRF) to 
ensure the coverage of a longer range. It is developed as is 
TDWR, to have high resolution and sensitivity [51] . Therefore 
some NEXRAD radars with a suitable windshear detection 
algorithm will be used in the airport terminal areas on an 
interim basis until the TDWR is deployed. 

In the detection of hazardous windshear conditions, 
reliable algorithms should be available to process the 
weather data obtained by a high precision Doppler radar. 
Considering a typical microburst characteristic, at least two 
methods are suggested. One is to compute the windshear hazard 
index which is represented in terms of spatial gradient of 
windspeed [11]. The other is to recognize the "S" curve 
characteristic associated with microbursts [12]. This "S" 
curve, showing mean windspeed versus range, develops as a 
strong downdraft induces an outburst of damaging winds on or 
near the ground. Usually this kind of information is obtained 
by estimating weather spectral moments through the DFT or the 
pulse pair method. Based upon these same basic concepts of 
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hazard detection, new techniques are developed and explained 
in this dissertation. 

Of course, the reliability of any detection algorithm 
depends on the quality of original data collected by a 
weather radar. Therefore the system phase stability problem 
is one of the important issues in developing windshear 
detection radars. Although the airborne Doppler weather radar 
is considered to have good potential for providing the pilot 
information to help avoid the hazardous windshear conditions, 
it may be susceptible to phase noise. The period over which 
the radar return can be considered stationary limits the data 
analysis window in an airborne weather radar as compared to a 
ground based Doppler weather radar. Also since the radar 
platform is moving there may be less advance observation time 
as the aircraft approaches the area of windshear hazards. 
These limited data records may yield poor estimates of 
autocorrelation or power spectrum parameters which will be 
further degraded by system phase noise. Phase noise may be 
compounded in an airborne weather radar which is operating 
under look— down conditions . Since clutter power is typically 
much stronger than weather return power, increased clutter 
spectrum width due to any system phase instability will 
obscure the weak weather return signal making the clutter 
rejection inherently more difficult. 
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Windspeed Gradient and F-Far.tnr 

The phenomenon of windshear is associated with air 
turbulence and can generally be characterized by large 
spatial gradients of wind velocity (windspeed gradient) . It 
seems feasible that the existence of a windshear condition 
may therefore be established on the basis of the measured 
wind velocity magnitude change within an airspace volume, 
given adequate spatial resolution of the measurements. With a 
typical pulse Doppler weather radar employed to remotely 
sense windspeed, the distribution of measured Doppler 
frequencies within each range cell is related to the 
reflecting particle motion within that range cell. The mean 
Doppler frequency of the return is a measure of the average 
windspeed within a range cell, while the spread of these 
frequencies is an indication of the turbulence within the 
range cell. 

Bowles [11] defined a hazard index F which is 
represented as 



where W x is the rate of change of the horizontal component of 

windspeed, v is the relative speed of the aircraft, g is the 
acceleration due to gravity and W h is the vertical component 

of windspeed. This F— Factor was derived considering aircraft 
energy balance for flight in spatially and temporally varying 
windf ields . Positive values of F indicate a performance 
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decreasing situation while negative values mean the opposite. 
In [11], it is shown that the average F values of all 
investigated aircraft accidents exceeded 0.15 which may be 
considered as a threshold value. Although the vertical 
windspeed is needed in the computation of F value, the 
forward-looking airborne radar can only measure the radial 
windspeed along the flight path of the airplane. Therefore 
considering only the measurable term W x , the radial component 
of F, F r is given as 

w; _ v aw x 
R g " g ar 

where Aw x is the change in radial velocity between range bins 
and Ar is the distance between range bins. From this 
expression, it should be noted that F R represents just the 
normalized value of windspeed gradient. Therefore, the direct 
measurement of windspeed gradient with low PRF radar may be 
very useful in the detection of hazardous windshear 
conditions. 

Detection and False Alarm Probability 
A reliable windshear detection system should alert a 
pilot when hazardous conditions exist, but nuisance alarms 
should be avoided by selecting a suitable statistic for 
detection. One fundamental problem is how to select the 
proper threshold value. A hazard index F-Factor has been 
suggested as appropriate [11]. Although a decision should be 
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based on both the detection and the false alarm probabilities 
along with risk factors associated with each possible 
outcome, it is difficult to compute those values since the 
probability distribution of F values must be known. If 
P (F | H-l ) and P(F|H 0 ), the conditional probability distribution 
of F in the presence of a windshear hazard, H lf and with no 
hazard, Hg, respectively, are known, a Nyman— Pearson 
criterion [13], for example, may be applied to set a 
threshold value which maximizes the probability of detection 
P DE t while constraining the false alarm rate P FA . An optimum 

test, thresholding the likelihood ratio, 

P (F | Hi) /P (F | H 0 ) =A (F) > X 

would indicate a windshear hazard alarm where P DET and P FA are 
given by 


P DET -Jp <Af Hl ) dA, P PA -Jp <A|h 0 > dA . 

A A. 

This likelihood ratio test is illustrated in Figure 5. In 
order to validate this test, numerous observations of 
microbursts should be made to estimate a reliable P(F|H 1 ) . 

This in turn can be used in determining a receiver operating 
characteristic curve which will show the relationship between 
the decision threshold value, P DET an< i P FA . An improved hazard 

index estimator has been suggested by Britt [14] using a 
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Figure 5. A Likelihood Ratio Test of Two Conditional 
Probability Density Functions 
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weighted least square technique to estimate windspeed 
gradient. This should provide a more reliable conditional 
probability distribution of F values with very little 
additional computation requirement. The basic component in 
the hazard index is the estimate of the spatial gradient of 
windspeed . 

Problem Statement 

Phase stability is a very important factor in obtaining 
accurate and reliable information in a coherent radar system. 
Since Doppler weather radar applications may impose severe 
requirements on phase stability, it is necessary to analyze 
the system phase instability problem which can be caused by a 
variety of sources, e.g., radar oscillators, complex sampling 
procedures, filtering and amplification, etc. [53] , [54] , [55] . 
Radar oscillator instability in a coherent radar has been 
previously analyzed, e.g, [15], [16], however, this work has 
been concerned primarily with a target detection radar 
system, not a weather radar system, and generally no analysis 
of the entire system phase noise effect is included. They 
have just considered one major phase noise source, i.e., 
STALO in relation with the performance of a single target 
detection radar. Therefore, in Chapter II of this 
dissertation, a general approach is developed to investigate 
the whole system phase noise effect. Two very popular weather 
signal processing algorithms, i.e., the DFT method and the 
pulse pair method, are analyzed considering phase noise. 
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Finally a benchmark laboratory test is proposed for the 
measurement of whole system phase noise in a pulse Doppler 
weather radar. 

The commonly used pulse pair method is quite attractive 
when processing an enormous amount of weather radar data in 
real time since it is considered the fastest algorithm 
available. The original concept of the pulse pair algorithm 
was developed by Miller and Rochwarger [8]. Since then, it 
has been studied and analyzed by many researchers. The pulse 
pair method was derived and has been evaluated most often 
under the assumptions that the weather spectrum is symmetric 
and relatively narrow. With the turbulent situations 
associated with windshear, these assumptions may not be 
valid. Some observed weather spectra [17] show that nearly 
25% are seriously skewed and can not be considered to be 
symmetric. This means that the pulse pair method may need 
reevaluation considering the skewness effect. Chapter III 
analyzes this effect using the skewed Gaussian spectrum 
mode 1 . 

The poly-pulse pair method was originally suggested as a 
way of enhancing the accuracy of spectrum moment estimation, 
but this method may be also useful in reducing the bias 
errors of a skewed spectrum. Based on the similar concept, a 
new modified pulse pair mean estimator is developed in 
Chapter III where it shows an improvement over a conventional 
method by reducing the bias errors. In the symmetric 
spectrum, the mean and the mode are same. However, in the 
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case of a skewed spectrum, it may be questionable that the 
mean is a more representative value than the mode. Hence the 
difference between the mode and the true mean due to skewness 
effect is also presented in Chapter III. 

Chapter IV presents a new method in the measurement of 
windspeed gradient. Since the windshear hazard index is 
computed from windspeed gradient, it may be considered the 
most important information in a windshear detection radar. 
Usually this information is obtained through the measurement 
of mean windspeed in each range cell. Because of the very 
strong wind involved with microbursts, a high PRF (pulse 
repitition frequency) radar is preferred to increase the 
Doppler spectrum range and thus the maximum unambiguous 
windspeed. This, however, results in an increased ambiguity 
in the range measurement. To overcome this difficulty, a new 
method is suggested for direct measurement of the windspeed 
gradient using a low PRF radar. 

Another challenging problem arises in dry microburst 
cases where a very low s ignal-to-clutter ratio may cause 
meaningless results without efficient clutter filtering. In 
some situations where clutter and weather spectra are not 
well separated, it is difficult to remove the clutter without 
deteriorating the weather signal. Even though efficient 
clutter filtering can be done using algorithms developed 
recently (12], [18], [19], the additional computation 
requirements may be prohibitive in real time weather radar 
signal processing. A new approach described in Chapter V does 
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not need any preliminary processing such as clutter 
filtering. It applies a modified Prony method to estimate 
peak points of weather return spectra. These values may be 
used to recognize typical "S" curve characteristic associated 
with a windshear phenomenon. Using a NASA simulation model, 
some validating results are also shown in Chapter V. Overall 
conclusions and some recommendations for future work are 
presented in Chapter VI. 



CHAPTER II 


PHASE NOISE EFFECTS ON WEATHER SPECTRUM 
MOMENT ESTIMATION 

Introduction 

In the weather radar return the superimposed scattering 
of incident electromagnetic energy from many randomly 
distributed particles causes an unavoidable phase jitter 
effect which generally contributes to the spectrum width of 
the weather radar return signal. With a coherent pulse 
Doppler radar, this creates some ambiguity in determining a 
representative windspeed condition within a particular range 
cell. Any radar system phase instabilities may also 
contribute to this ambiguity and affect the Doppler frequency 
resolution as well as the dynamic range capability. There are 
many potential sources of system phase instability (phase 
noise or phase jitter) but a potential primary source is the 
stable local oscillator (commonly STALO) which provides 
transmitted and reference carriers within the radar. If the 
range time is adequate to decorrelate the STALO oscillator 
phase, any phase variation with time (jitter) may be adequate 
to contribute error in return Doppler spectral estimates. 
Furthermore, this may not be apparent when observing the 
weather return Doppler spectrum because of its inherent 
spread. 
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This chapter provides an analysis of the effect of pulse 
to pulse system phase jitter uncertainty on the pulse pair 
[29] and DFT based spectral moment estimates. More 
importantly a general analysis approach is developed which 
will allow for a quantitative assessment of the effect of the 
entire radar system phase noise. The approach can be used 
with an analytical model description of phase noise or with 
an actual measured system phase noise power spectrum. 
Intrapulse phase uncertainty is not considered here. 

An analytical development of the Doppler spectrum using 
a Gaussian approximation for the phase jitter spectrum is 
presented first. Even though a power-law model is more 
appropriate in representing experimental results of a STALO 
phase noise spectrum [20], a Gaussian model is considered 
here for two reasons. First, it is reasonable to fit a 
Gaussian distribution to a major portion of a power-law phase 
noise spectrum particularly in the tails of the spectrum. 
Secondly, a Gaussian assumption provides mathematical 
tractability . It is recognized here that spectrum moment 
estimation error may depend upon the specific model of phase 
jitter noise. It is also true that the optimum model is not 
known. For example, considering an oscillator alone, even 
though many related papers have been published, none of them 
exactly describes all instability factors [2 1 ] , [ 22 ] , [ 23 ] . 
Therefore, the approach developed here is generalized to 
enable analysis of all types of phase jitter effects through 
numerical computation of the autocorrelation function from a 
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given measured or otherwise specified noise spectrum. This is 
particularly important because the analysis procedure can be 
used for any set of phase jitter power spectrum measurements 
that might be available or for any analytical spectrum 
function defined in terms of a set of values. The spectral 
mean and width as estimated by the pulse pair method and DFT 
method are analyzed analytically considering the effect of 
phase jitter. For the Gaussian model, the effect of phase 
jitter is computed for a particular phase jitter power level 
while allowing variation in both the weather spectrum width 
and the phase jitter spectrum width. 

Gaussian Phase Noise 

Any phase jitter occurring at the transmitter stage will 
spread the original sinusoidal signal spectrum. The stable 
local oscillator output is modelled as 


V(t)=V 0 COS(o> c t+<t>(t)) 

where <l>(t) is a phase jitter and the carrier is CO c =27Cf c . If 
one assumes that O(t) is a normal stationary process with 
zero mean and variance a 2 = E[<I> 2 (t)], the spectral density of 
the transmitted signal is given by [24] 

2 

s.(0=^ ! (S(f-f ) + s (f-y + £-l [s (tjsk (f )] f ) 

n=2 c 


(2 . 1 ) 
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where the first term inside the braces of (2.1) relates to 
the sinusoidal energy at f c , the second term is the phase 
noise spectrum located at f c and the last term includes n-1 

convolutions of the phase noise spectrum with itself, 
followed by a translation around the carrier frequency f . 

As the mean intensity of phase noise becomes very small, the 
spectral density of phase jitter can be approximated by a 
Gaussian spectrum [25], [26]; i.e. 

S /f > w 0 e*p{lJ) 

c 


where the relationship between phase jitter width Aco c 2 and the 
total phase jitter power o 2 is represented by 


oo 



With Gaussian phase noise , the spectral density of the 
transmitted signal from (2.1) can be obtained analytically as 


S v (f)= 





exp[- 


-(CO-CD f 

— r 1 ) 

n Aar 

c 
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If (J is small (a< 0.5 radians), which should be true in a high 
quality oscillator, S v (f) can be approximately given by 


S (f)=^e '° 2 
v 2 


{ 5(f-f) 


2 ^ 

Aco 


a 2 exp[ • 


-(CO-CO J 2 


Aco 


]} 


( 2 . 2 ) 


It can be noted that the theoretically exact description for 
the spectral density is [27] 


vi 


s v (f)= -y J ex P(-y 1 2 ^)) cos( 0 ) c x) e' J<OT dx 


where 



E{[ <Ky^H(yf) 
(w c x ) 2 


This theoretical expression would be useful if an 
approximation for the true variance I 2 (X) could be determined 
for all X, but it is very difficult to obtain this 
experimentally except for a few discrete values of X due to 
the requirement of very long observation time. Therefore, the 
approximation in (2.2) provides a realistic means of 
analyzing the effects of phase jitter in the Doppler radar. 
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The transmitted radar signal is scattered by many 
particles and the signal frequency will be shifted according 
to the particles' velocity which is the well known Doppler 
effect. Using the description given in (2.2), the scattered 
return signal can be approximated as 


S v (f) = 


V , 

= — ® e^C { - z 
2 ^ 


-exp (• 


W 


(f - f c~ f d ) 

2W 2 


-) to 2 - 


yj 27C(Af 2 /2+W 2 ) 


• exp ( 


(f-fc-f d > 2 
2 (Af 2 /2+W 2 ) 


) } 


(2.3) 


where W is the spectrum width with zero transmitter phase 
noise and C is a constant determined by the reflectivity of 
target particles and some other factors (see for example 
[28]) . If the noise is small (<J small), the second term in 

(2.3) can be ignored and the typical Doppler return signal 
will have a purely Gaussian spectral shape. With this 
assumption the autocorrelation function of the return signal 
after demodulation will be 

2 

V 

RW=YCKo,w,i)e W +N6 t0 (2 4) 


where 
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-y(a,w,T)= 



-. 22 , 

-2n w t z , 

e +a z e 


„ 2 . 
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The additive white noise N Sr, o is included in the signal model 
to represent the usual background noises [29]. The 
autocorrelation function shown in (2.4) is obtained under the 
assumption that no phase noise occurs at the receiver in the 
process of demodulation. However, in the interest of a more 
complete understanding of pulse pair estimation, the effect 
of phase jitter occurring at the receiver should be included. 

In the pulse-pair method [30], a two-point estimate of 
the complex autocorrelation is based upon processing coherent 
returns from pairs of transmitted pulses. Assuming that the 
pulse pairs are separated by T s and repeated with period T as 

illustrated in Figure 6, the complex time autocorrelation 
function of the complex received wave form z(t) at the lag 
value T s is defined as 

R(Ts)=E { Z*(iT) Z(iT+T s ) ) . 



Figure 6. Pulse Pairs for the Estimate of Autocorrelation 
Function 
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Considering the phase deviations in the process of 
heterodyne demodulation designated by O iT and O iT+Ts , the 

autocorrelation function is 

R(T)= R(T) Ele'S*'™'! (2 . 5) 


where, from (2.4), 


R(T S )» S 7(a,w,T s ) e J< ° dT ‘ 



(2.6) 


If one assumes that and <!> 2 =<I> iT+Ts and are, for example, 

normal random variables [20], then (2.5) can be rewritten as 


R (T )= R(T ) ) 1 JJe J ‘"'c* 1 exp[-«t^+^-2p<T s )<|> i< ti. 2 ) 


A24>’(l-P 2 ),] d<t> d* 2 


where 


p(T s )= E^] / <j> 2 and ^ . 


The complex autocorrelation function can be further reduced 


to [31] 
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R(T)= R(T) exp[-W 2 -<t> x 4 > 2 )] 

oo 

= R(T s ) exp[J 1 -COSCOT ) df] 

-oo 


(2.7) 


Similarly S ( j,(f) can be assumed to have a Gaussian noise 

spectrum as described in the previous section since the same 
oscillator signal is fed into the receiver for heterodyning. 
With this assumption (2.7) will become 

R(T ) =R(T ) mTH N5 t Q (2.8) 

-AcoV 

'PCo.T )= exp[-a 2 { l-exp( ^ 5 )}] 

where R(T S ) is given in (2.6) and <T is the r.m.s. phase noise. 

As seen from (2.8), this modified autocorrelation at lag 
value T s is explicitly a function of the phase spectrum width 
Aco c and the total phase jitter power a 2 with the additive 
white background noise also included. 

General Approach by Numerical Computation 
Using the system approach, the general complex 
autocorrelation function of the weather return signal in the 
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presence of phase noise can be obtained without full 
analytical representation of the phase noise spectrum. Again 
representing the weather return Doppler spectrum as Gaussian, 
modulation and demodulation procedures also can be 
represented by a simple mathematical form considering phase 
jitter as can be seen from the Figure 7. Rj and R w are 
autocorrelation functions of phase jitter. Rj can represent 
the effect of phase instability in the radar system. R w 
represents phase jitter due to particle scattering. Following 
the same procedure used in deriving (2.8) [26], these 
autocorrelation functions can be written as 


Rj=E[e e ] =exp(-{Js^(f) df- R^CT )}) 



C : attenuation constant 


where 

oo 

f j2*fT f 

v^Jv 06 * 

-oo 


(2.9) 


Since the exponential term in Rj represents the difference 

between phase noise power and the autocorrelation of phase 
noise at T s , it can be shown that Rj goes to 1 . 0 (no phase 

noise effect) as the phase noise spectrum becomes narrower. 
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Considering statistical independence between these terms, the 
output complex autocorrelation function will be just the 
mu ltiplication of all individual autocorrelation functions 
which is shown to be 


• f -2 ji 2 wV jco/T 

R(T )=Sexp(-2{ S (f)df-R (T)})e e + N5 t n 

JO 0 s T,,0 


As seen from (2.9), a closed form expression for R^CTg) is 

generally not possible unless the phase noise spectrum is a 
known functional form, e.g., Gaussian. Therefore, some 
numerical method may be needed to determine Rj _ An adaptive 

quadrature algorithm [32] is used here in the computation of 
these terms. The next section analyzes the effect of these 
modified autocorrelations on the quality of the pulse pair 
and DFT estimates. 


Pulse Pair Estimation Errors 

Modern radar oscillators are very stable with small 
phase noise power so that the effect of phase jitter in a 
return Doppler spectrum may not be apparent, even if it does 
have an effect on the pulse pair estimates. Figures 8, 9 and 
10 illustrate this situation. In Figure 8 simulated return 
Doppler spectrum is shown with and without phase jitter. A 
zero mean normal signal return spectrum with unity total 
power and a width corresponding to 15% of the processing 
bandwidth is depicted along with the resulting spectrum 



Doppler Frequency (Hz) 


with phase jitter without phase jitter 


Figure 8. Example of Two Simulated Return Doppler 
Spectra Considering Phase Jitter Effect 
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with phase jitter ; without phase jitter 


Figure 9. Error Standard Deviation of the Spectrum 
Mean Pulse-Pair Estimate for the Example in Figure 
8 




with phase jitter 


without phase jitter 


Figure 10. R.M.S. Error of the Spectrum Width 
Pulse Pair Estimate for the Example in Figure 
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contaminated with statistically independent phase noise 
having a total power level at 9% of the signal power and a 
width which is 28% of the processing bandwidth. This phase 
noise width is less than that of an inverse f type 
distribution with the same total phase noise power. There is 
little apparent difference between the two spectra depicted 
in Figure 8. Analyzing the standard deviation of the error 
of the pulse pair mean estimate for the situation depicted in 
Figure 8 (wT s =0.30), however, shows in Figure 9 that the 
phase jitter causes more than a 15% increase. For other 
Doppler spectrum widths, this same phase jitter error can 
cause up to 50% increase in the error standard deviation 
without an apparent change in the spectrum. From Figure 10, 
the rms error of the pulse pair width estimate also can be 
significantly affected by this phase jitter, causing 
approximately 15% error for the Doppler return spectrum width 
considered. These results suggest that it may be useful to 
analyze and quantify the phase jitter effect on pulse pair 
estimation errors even when the effect on the Doppler 
spectrum is not apparent. 

Consider estimation of the mean frequency with pulse 
pair processing [29]. The estimated return spectrum mean is 
defined by 


( 2 . 11 ) 
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where 


R(T S ) 



i=0 


(iT) Z ( iT+T s ) 


As seen from (2.9) and (2.11), the estimate of the mean 
frequency will not be dependent upon phase jitter because it 
does not affect arg{R(Ts) }, i.e., only the magnitude scaling 
is changed by phase jitter. Thus, according to this 
analysis, no bias will be introduced into the mean frequency 
estimate because of phase jitter. 

The variance of the mean frequency estimate as 
previously published [30] is 


VAR(f d ) = [87t 2 Tj3 2 (T s ) ] _1 {m' 2 [1-P 2 (T s ) ] .^J3 2 (mT) (M- | m | ) 


m=- (M-l) 


+£5*gu-l»«T.>8 ; . V0 ♦ 


,) 

M T-T,,0 J 1 


( 2 . 12 ) 


This expression can be rewritten considering phase jitter by 
replacing (3(T S ) with (3’ (T s ) where 

P(T s )=7(a,w,T s ) XK(a,T ) 


for the Gaussian spectral model and 
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2ji z w'T 

P(T )= exp(-2{ IS (f) df-R (T)})e' 
s J 4> 4» s 

- oo 

for the arbitrary spectral model. Here M is the number of 
sample pairs. 

Analyzing the variance expression in (2.12), if the 
total phase jitter power a 2 is held constant, one can 

investigate the effect of varying the spectral distribution 
of the phase jitter. The parameter Afc determines the shape 
of this spectral distribution. As Afc increases, the phase 
jitter spectrum is broadened. As can be seen from Figure 11, 

which is a plot of the mean frequency estimate standard 
deviation versus normalized Afc, when the phase jitter 

spectrum is broadened, an increase in the variance of the 
mean frequency estimate can be expected. Also from Figure 
11, as the weather spectrum width W increases, this effect on 
the mean frequency estimate is more pronounced. Figure 12 is 
included to show the effect of phase noise through numerical 
computation from experimentally obtained phase noise data 
[34] . From Figure 12 it can be said that the Gaussian jitter 
spectrum model yields a maximum error bound in the estimation 
of mean frequency. 

Considering the weather spectrum width estimate using 
the pulse pair method, the width estimation is independent of 
the mean frequency [9] . If the width is sufficiently smaller 








no phase noise 

Gaussian model 


numerical example 


Figure 12. Comparison of Error Standard Deviation of the 
Spectrum Mean Pulse Pair Estimate 
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than the Nyquist interval, one form of the estimator which 
is not dependent upon spectrum shape is given by [9] 



(2 . 13) 


where 



k=0 


and Sis an arbitrarily small number. With the introduction of 
phase jitter, this estimator becomes 


W= \ 


/2 


JIT* 


R(T S ) 

A 

s 


A V 

when R(T S ) 


<S and 


when R( T S ) >S 


(2.14) 


where since |R(T S )| from (2.7) replaces |R(T S )| in (2.13), the 

estimate of spectrum width is biased. Previous results for 
the asymptotically unbiased variance can now be extended to 
yield the variance of the width estimate as 
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VAR(W) = [32M7t 4 (WT s ) 2 13 2 (T s )T;]" 1 {2 • [1-(1+8 T _ T Q )P* (T, 


+5 t _t 0 P 4 (T s ) ]~ + [ 1+ ( 1+5 t _t„ 0 ) P 2 (T s> l- 2 

3 


T-T , O' - ' ~ s ' J g 

+p 2 (T s )^{2p 2 (mT) + 

m=- (M-l) 


p 2 (mT) p 2 (T s ) 


+p(mT+T s (l-8 T _ T >0 ))-4p(mT+T s ) .p (mT ) f 1 (T s ) 


m 1 
<j — » > 


where P(T S ) given in a previously derived result [9] is 
replaced by P' (T 3 ) as was done in (2.12). The variance of the 
width estimate may be increased by phase jitter, but does not 
completely describe the goodness of the width estimator, 
since phase noise will bias the width estimate of spectrum as 
seen from (2.14) . Thus, in this situation the mean square 
error defined by 

E { (0-0) 2 }= VAR (0 ) + B 2 


where B is the bias of a biased estimation 0 of the parameter 
0, will be a more representative measure of the effect of 

phase jitter on the pulse pair estimation of spectrum width. 
The bias term and the root mean square error of the width 
estimate are plotted in Figures 13 and 14 respectively. 
Figure 14 shows that the variance of width estimate is also 
increased due to phase jitter. In each case the normalized 
weather spectral width has been varied to illustrate its 
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0.40 


0.05 


Figure 13. Bias of the Spectrum Width Pulse Pair 
Estimate Considering STALO Phase Jitter 
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effect upon the estimation errors. As is seen in Figure 13, 
the bias in the width estimation is more sensitive to 
increasing phase jitter when the true normalized spectral 
width is small. This variation in sensitivity is less 
evident when considering the r.m.s. error in the width 
estimate in Figure 14. It should be noted that in each of 
Figures 11, 13 and 14 the curve for Af c =0 (no phase jitter) 

has been previously published [30] . 

DFT Estimation Errors 

Consider the problem of calculating the mean frequency 
and the spectrum width from the DFT-derived spectral density. 
The most straightforward DFT spectrum estimates are simply 
the mean and width estimate of the derived spectral density 
S(fi) [7], i.e., 

f- -i 

i 

These, of course, are parameter estimates of the composite 
signal-plus-noise spectrum. As it is well-known, these are 
not satisfactory estimates of signal parameters except for 
large signal-to-noise ratios. Therefore, some technique of 
removing the noise influence should be implemented. One 
method of noise suppression is to simply subtract the noise 


(f.-f) S(f 1 ) 


and -- 2 - 


w = 


X s < f i> 

i 
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spectral density N(f^) from the derived spectral density and 
calculate the estimates of the resulting spectrum, i.e, 

X f i [S(f i )-N(f i> ] X (f i -f) 2 [s(f i )_N(f i> J 

1 '2 i 

f- — — ; w 

XtS(f i )-N(f.) ] 21 [S < f i> “ N < f i> ] 

i i 

These estimates are unbiased by noise even for low signal-to- 
noise ratios. This method has been chosen here as the basis 
for evaluating the phase noise effect since the variance of 
both these estimates was derived earlier by Berger and 
Groginsky [33] . Another variation [7] to eliminate biases due 
to aliasing is useful for symmetric spectra, but it is not 
adopted in this analysis because no variance expression is 
available. Simulation results do, however, show similar 
estimation quality when compared to the method chosen here 
for further analysis. 

Using previously derived variance expressions [33] and 
replacing the summations in the estimates with an integral 
over the Nyquist interval in each case, the mean frequency 
estimate variance can be expressed as 



1 

2T 


1 

2T, 
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var (f ) = 

2 

MT 


s 


f 2 2 

[T f S (f+f ) 

sj n m 


1 



df + 2 ( 


N 2 f 2 

>t f s 
s S J n 


(f+f ) 

m 



df 


+ 


1 N 
12 < S ) 


2 

] 


and the width estimate variance is 
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) } 


where 

f : the true mean frequency, 

S n (f): the normalized power spectrum, i.e, 

S n (f)=S (f L ) /S where S=]Ts(f.), 

W: the true spectrum width and 
N: a constant white noise power . 

For a narrow Gaussian spectrum without phase noise, the 
variances can then be simplified to 
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Unless the spectrum width is small compared to the Nyquist 
interval then zero mean assumption is needed to validate 
results. With the zero mean assumption, this analysis of the 
phase noise effect can validate the simplified derivations 
for other than very narrow spectra. 

Considering Gaussian phase noise. The normalized 
spectrum S n (f) becomes 


S n (f)= e 


-2G 


2 V* 


i— 2 Af 2 7 

J 2n (w +k -) k! 

2 


exp (- 


2 Af 

2 (w +k — ) 

2 


(2 . 15] 


where w is the width of the Gaussian spectrum without phase 
noise. Since the phase noise power is very small ( <J< 0.5) in 

quality oscillators, the normalized spectrum can be 
approximated very accurately by the first four terms, i.e.. 
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Using this spectrum expression, the variances were rederived 
and the results are given in Appendix A. From these 
expressions, one can investigate the effect of varying the 
spectral distribution of the phase jitter (parameter Af c ) and 

also the effect of varying the phase noise power (parameter 
<* ) * Figure 15 shows that the standard deviation of the mean 
estimate is increased as Af c T s increases for the given phase 

noise power. From Figure 16, it also can be seen that 
broadening the phase noise spectrum increases the standard 
deviation of the width estimate. 

Although a Gaussian phase noise assumption may be 
reasonable, and with this assumption a mathematically 
tractable analytic analysis is possible, the effect of phase 
noise may actually be dependent upon the specific noise model 
selected. For this reason. 


as in the previous section. 
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another approach was developed to analyze the effect of any 
given phase noise spectrum S<j)(f) which may be defined by 

interpolated experimental data or which may be some other 

kind of appropriate phase noise model. For the specific phase 
noise spectrum, < f ) , the procedure to obtain the variances 

of Gaussian Doppler spectrum with phase noise is described in 
Figure 17 using block diagrams. With this approach, published 
interpolated phase noise data [34] were processed to 
investigate DFT based spectral parameter estimation errors. 
As seen from Figure 18, the phase noise effect is almost 
negligible as might be expected, since the given phase noise 
spectrum is extremely narrow and contains very little power. 

Limitation of Clutter Cancellation Capability 
Ground clutter returns in an airborne Doppler weather 
radar are typically dominated by a very strong spectral power 
around the Doppler frequency corresponding to the aircraft 
ground speed. This is treated as the zero frequency reference 
in the discussions presented here. This ground clutter return 
is assumed to be from stationary structures on the ground and 
the ground itself. Additional clutter may also be present as 
return from moving objects (vehicles, etc.) on the ground. 
Without proper clutter cancellation, the spectrum moment 
estimates associated with the weather return may be so 
seriously biased that they do not provide reliable 
information for windshear detection [35] . It is considered 
that a 10 dB s ignal-to-clutter ratio (SCR) is needed for 
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accurate mean velocity estimation via a pulse pair processor 
while the spectrum width estimation may require 15 dB SCR 
[36] . This means that at least a 10 dB SCR should be achieved 
by filtering clutter power. Clutter power located around zero 
frequency can be removed very effectively using a high pass 
filter which attenuates low frequency powers to the desired 
level in a stopband region, but the width of stopband should 
be carefully selected to avoid elimination of low frequency 
weather signals. When the clutter spectrum is narrow, clutter 
cancellation can be done very successfully without degrading 
the weather spectrum. However, phase noise in a radar system 
will broaden the clutter spectrum making some portions of 
power spill into a pass band region limiting clutter 
cancellation capability. This may be a problem particularly 
in the case of a dry microburst situation where the weather 
signal is much weaker than the clutter return. 

If it is assumed that the clutter and weather return 
spectra are both Gaussian, the return spectrum can be modeled 
as the sum of two normal functions 

C f 2 1 

S(f) = -— — exp ( ) + -—Z — exp ( — ) 

V2JIW C 2w c V27 IW s 2W 2 (2.16) 

where C = 10 , f^ is the mean frequency of a weather 

return spectrum, the clutter is assumed to be zero mean, and 
W c and W 3 are the spectrum width of clutter and weather return 
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respectively. Considering the phase noise effect and using 
Equation (2.15), Equation (2.16) is rewritten as 


S (f) = Ce 


k=0 


2 V k 


/lit (wf+k 
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-exp (■ 


c. 1/2 


O / TYT _L 1,. 
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+ e-^£- 
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k -°/^ <w| + k %> 1/2 
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2 

2<W s +k— -) 


-) 


(2.17) 


Here, it should be noted that a Gaussian phase noise spectrum 
has been assumed. Since G is usually very small as described 

before, Equation (2.17) can be approximated using the first 
four terms of the infinite summation. Considering the use of 
an ideal high pass filter having the stopband width of 2B and 
an attenuation value of ATT in dB is assumed, the SCR of the 
filtered return spectrum is represented by 

SCR = 101og 10 <!f) 

where CR and SR can be expressed as 

C R = Ce' 202 [Q (-f-) +2o 2 Q (-f-) +2c 4 Q <-|-) -4<J 6 Q ) ] 

w c o 

+ 10 -ATT/10 Ce _2o 2 [ (1 _ Q( _|_ ) } + 2 q 2 ( 1-Q (-£-) ) 

w c W 1 

+2a 4 (i-Q(-|-) )4-|a 6 (i-Q(-^-) ) ] 

W 2 3 W 3 
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The resulting SCR figures after filtering the return spectrum 
are plotted for various values of f d by changing phase noise 
power a 2 . It has been shown elsewhere [35] that the SCR of a 
filtered spectrum is decreased if the weather spectrum is not 
well separated from the stopband region. The emphasis here is 
placed upon the phase noise effect. It should be noted that 
the results are obtained assuming use of an ideal high pass 
filter and no aliasing effect in terms of the spectral 
representation. Figures 19 through 22 are obtained by 



SCR (dB) 



phase noise power 


Figure 19. Phase Noise Effects on Clutter Filtering 
of 50 dB Stopband Attenuation for a Doppler Weather 
Return with 0 dB SCR 
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stopband width=150 (Hz), ATT=50 dB 



Figure 20. Phase Noise Effects on Clutter Filtering 
of 50 dB Stopband Attenuation for a Doppler Weather 
Return with -30 dB SCR 


SCR (dB) 





SCR (dB) 
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changing the attenuation value with the fixed stop band width 
of 150 Hz. All plots show that the SCR of a filtered spectrum 
is decreased as the phase noise power increases. The 
undesirable effect of clutter filtering resulting in some 
elimination of the weather signal power is also noticed. As 
expected. Figure 19 shows that phase noise is not a serious 
problem if the weather signal is strong or comparable to the 
clutter return. However, as can be seen from Figures 20 and 
21, any small phase noise may seriously limit the clutter 
cancellation capability when the clutter power is much 
stronger than the weather signal. Even with the ideal high 
pass filter with 70 dB attenuation in the stopband, more than 
10 dB SCR cannot be achieved if phase noise power exceeds 
0.05, as shown in Figure 22. 

ExQPPSed Measurement of Radar System Phase Noise 
It is desirable to determine the effect of phase jitter 
noise on the error in the pulse pair and DFT parameter 
estimates to understand the expected radar system 
performance. This section describes a procedure for the bench 
test of the radar to yield an estimate of the effective phase 
jitter power spectrum. For this purpose, the transmitted 
signal with some delay time is to be fed directly into the 
receiver. I,Q data would then be obtained to analyze phase 
jitter in the absence of a return Doppler signal. A return 
signal is not included in the procedure since just the pure 
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phase noise term is of interest. Figure 23 illustrates the 
setup . 

It is anticipated that the stable local oscillator 
(STALO) will be a major potential contributor to phase 
jitter, even though it may be very small. The following 
analysis is based upon modelling the STALO phase as time 
varying. At the receiver illustrated in the Figure 23, the 
demodulated noise signal is 


V Q (t)=I (t) +jQ(t) = 


e 


-0(t-t d ) ) 


where t d is the delay line time (delay time is less than 25 
|lsec in the case of 7.5 km maximum detection range) . Here, 

STALO oscillator phase jitter can affect the modulated signal 
at the transmitter and after a delay t d can affect the 

demodulated signal. The interest here is in measurement of 
the effects of pulse to pulse phase jitter as will be 
observed at the receiver. 

The autocorrelation of the receiver output Rj(T s ,t d ) will 
be a function of both the transmitter/receiver delay time t d 
and the pulse pair separation T s . It can be expressed as 
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R j (T s' t d )=E[V o (t)V o (t+T s ) ] 

j (0 (t ) -0 <t-t d ) ) - j <0 (t + T g > -0, (t+T.-t d ) ) 

=E[e e 1 1 ] 

j0.(t) -j0 <t +-T S ) -j0 <t-t d ) j0.(t + T s -t d ) 

=E [e e e e ] 


(2.18) 


where <t>i (t) represents the phase variation with time in the 
radar oscillator. With a PRF on the order of 3-4 kHz, the 
pulse separation is perhaps two orders of magnitude greater 
than t d . It is expected that phase variations within the radar 

would be very slow such that the intra— pulse variation can be 
considered insignificant, but the delay time t d may be long 

enough to experience a measurable phase variation. Thus if 
t d <l |lsec the radar phase at the transmitter and at the 

receiver can be considered identical for each pulse and 
jitter present at the transmitter is exactly cancelled at the 
receiver, i.e., Rj=1.0. In the case of much longer t d delays, 
<{>l(t) and <t>i (t— t d ) can be considered uncorrelated and the 
expression in (2.18) can be reduced to 


Rj <T ,oo)= Ri (T s ) =E[e 


j '-V 


j0 1 (t) -j0 <t+T s ) -j0 1 (t-t d ) j0 1 (t-t d + T s ) 

e ] E [e e ] 


j0 (t) -j0 (t + T s ) 2 

= {E[e e ]} 


:-2 ( J £ 


= exp(-2{ SJf) df-R (T ) } ) 


(2 .19) 



66 


using (2.10) . Equation (2.19) represents the worst case that 
can be assumed since no phase cancellation occurs. With 
realistic delays, t d , what actually will occur is probably 
Rj (T s ) <Rj (T s , t d ) <1 . 0, since Rj(T a ,t d ) would generally decrease 
as t d increases . 

The analytical result suggests that a practical 
measurement of Rj(nT s ,t d ) may be possible using sampled I,Q 

data if the number of samples is very large , i.e., 

M-l 

i v 

(nT )= — / ,V ( iT) V (iT+nT ) 

J s 0 o s 

i=0 


where nT s is the delay between pairs of pulses separated by n 
pulse repetition periods. Characterizing R 1j (nT s ,t d ) versus n 
for specific values of t d will then delineate the phase jitter 
pcwer spectrum due to the radar instabilities. The resulting 
measured complex autocorrelation could then be used in (2.12) 
and (2.13) to correctly specify the parameter estimation 
error for the actual system used. 

Of course, it is expected that the measurement results 
will represent Rj as a composite of all system phase noises. 

Summary 

This chapter has developed a procedure for specifying 
the phase stability of a Doppler radar system through 
definition of the relationship between pulse to pulse phase 
error (phase jitter) and weather spectrum parameter 
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estimates. Both the pulse pair and DFT based estimators have 
been considered. The effect of phase jitter on post clutter 
rejection filtering signal to clutter ratio has also been 
analyzed. Finally a procedure for measuring the phase 
stability of a radar system has been proposed. 



CHAPTER III 


SPECTRUM SKEWNESS EFFECTS 
ON PULSE PAIR ESTIMATION 


Introduction 

The pulse pair estimator calculates the first two 
moments of the Doppler spectrum from estimates of the complex 
autocorrelation function at lag T s Goodness of this 

estimator is typically determined by examination of the bias 
and the variance of the moment estimates. To analyze the bias 
in the pulse pair estimates, consider the process 
autocorrelation function R(T S ) expressed in terms of the true 
mean Doppler frequency f^ [28] : 


2T S 

) 2 * f d T s f . i 2ItfT s 

R(T s )= e Is (f ) e df 

-l (3.1) 



where s'(f) is the zero-mean representation of the weather 
Doppler spectrum. Unbiasedness of pulse pair estimates is 
based on the assumption that a spectrum is symmetric or so 
narrow that the imaginary part of the integral in (3.1) can 
be considered as zero, i.e. [28], 
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l 



Js'(f) sin (27tfT s ) df = 0 


-l 



However, the weather return Doppler spectrum is often broad 
and not symmetric thus causing a bias. This bias effect is 
analyzed here using a skewed Gaussian spectrum model with 
various spectrum widths. 

Analysis of Bi as Errors 

In this analysis, a skewed spectrum will be modelled as 
piecewise Gaussian with appropriate normalization, given by 



when f<0 


when f>0 


(3.2) 


where the standard deviation ratio p=wi/w 2 defines the degree 
of skewness, g, i.e. [37], 



[ (p ‘ 2 + l) 


3 3 

2 2 2 

"(P +D ] . 


g= 
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This skew parameter varies proportionally to skew from g=0 
for no skew (p=l) to larger values, e.g., g=3.14 for a case 

which may be considered large skew (p=10) . Figure 24 shows 
the relationship between the parameter, p and the degree of 
skewness, g. For a narrow Gaussian spectrum with symmetry, 
i.e., w^=W 2 =w, the integral in the autocorrelation function 

(3.1) can be reduced to one simple term, exp {-2n^vi 2 T s 2 ) , but 

the results for the skewed spectrum model will include both a 
real term 

l 

2T s 

a= S (f ) cos (2JtfT ) df 

J n s 

1 

2T S 

and an imaginary term 


2 
f 

2 

2 

e ) ] sin (2lCfT ) df 

S 

Using these terms, the bias in the pulse pair mean and width 
estimates can be represented by 



222 222 
2 1 -2JI Wl T s W 2 -271 w 2 T s 

1 +p ^2 2w x C 





72 


mean bias= 


- — tan 


271T 


- f„ 


width bias= 


ll 7tT c 


In 


l 


VVa 2 +b 2 ; 


1/2 


-w 


(3.3) 


where true mean, f m and true width, W are described from 
Equation (3.2) as 



f = 

m 


fS (f) df= 

n 




(pw 2 -w x ) 


1 

2T_ 


w 


■J 


<f- 


f J 


1 

2T e 


1 2 2 2 

: f ) df= — — (w +pw ) -f 

i+p 1 2 m 


Estimate biases as given by Equation (3.3) are plotted as 
functions of the true width W and the skewness parameter g in 
Figures 25 and 26. In Figure 25 if there is no skew (g=0) tne 
pulse pair estimator is unbiased. As skew is increased there 
is a sharp increase in the bias. Once the skew parameter g>0, 
the percentage bias error is essentially independent of the 
specific value of skew but is strongly related to the 
spectrum width W. As seen, the bias error due to skewness is 
not negligible if the spectrum is broad. Figure 26 shows that 
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Figure 25. Mean Bias Error in Pulse Pair Estimates 
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Figure 26. Width Bias Error in Pulse Pair Estimates 
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a broad spectrum with a large degree of skewness can degrade 
the quality of pulse pair width estimates, but it does not 
seem to be as serious compared to the pulse pair mean 
estimate error in Figure 25. 

To get a more complete measure of the effect of skewness 
on the pulse pair mean estimate. Figure 27 compares the 
estimate r.m.s. error for the case of the skewed spectrum 
with that of a symmetric Gaussian spectrum having an 
equivalent width. As seen from Figure 27, the error caused by 
the skewness may seriously degrade the pulse pair estimation 
quality if the return Doppler spectrum width is 40% or more 
of the Nyquist bandwidth. 


Polv-pulse pair Method 

Since the pulse pair mean estimator bias is sensitive to 
skew, an alternative may be desirable. Therefore the poly- 
pulse pair method [38] is investigated here as a possible way 
of minimizing such errors. The pulse pair mean estimator 
algorithm uses the first lag of the complex autocorrelation 
estimate and is based on a linear approximation to the 
derivative of the phase function 
autocorrelation estimate, i.e., 


a = _l_ d0CT s ) | = _i_ 
271 dT c T s =0 - 2 71 


0(TJ 


of the 


comp lex 



R.M.S. Error of Skewed Spectrum 


R.M.S. Error (standard deviation) of 
Symmetric Spectrum with Equivalent width 


Figure 27. R.M.S. Error of Skewed Gaussian Spectrum 
for a Certain Degree of Skewness (g=1.99) 
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A 

where 0(T S ) is the phase function. There will be no 

approximation error for a symmetric spectrum, but a large 
error can occur in a skewed spectrum since 0(T S ) is no longer 

a linear function of T g . An alternative which may reduce 

these errors is to approximate 0(T S ) as a low order polynomial 

(greater than first order, n=l) , i.e., 

® tr - ) ,3.4, 

i=l , odd 

where 0(T S ) is an odd function of T s [8]. An odd function 
representation is needed since R(T 3 ) is the Fourier transform 

of a real valued spectrum. Then the mean estimate of a skewed 
spectrum will be 

- d ^s) | = A 

2 k dT s t s =o “ 271 a i ' 

To estimate a^ for a particular n>l, the complex 

autocorrelation function must be estimated for lags other 
than the first lag value T s . This can be accomplished with 

the poly-pulse-pair method [38] . Figure 28 shows that the 
mean bias error for a spectrum skew of g=1.99 (p=2) can be 

significantly reduced over a range of spectrum widths using 
the poly-pulse-pair method. To more completely evaluate the 
poly-pulse-pair method the variance of these estimates should 
be compared to the conventional pulse pair method (n=l) . 
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Considering a third order polynomial model in Equation (3.4), 
i.e . , 


ecr s ) = + a 3 T s 

the poly-pulse-pair mean estimate variance can be shown to be 


2 

A A T- 

A A 

T 4 


var' (f) = var(f) - E 

a a -a a 

3 1 3 1 

S , E 

( a -a 1 

2K 2 


47I 2 

[ 3 3 ) _ 


where var(f) is the variance of the conventional pulse pair 
method. The first term may be positive or negative and the 
second term will actually reduce the pulse pair estimate 
variance. In any case, since the pulse interval T s is 

generally very small, the higher order terms may be ignored 
to yield 

var ' (f ) = var (f ) • 

From these results, it appears that the poly-pulse-pair 
method can improve the quality of mean estimates by reducing 
bias errors in a skewed spectrum. 

Mode versus Mean Estimation 

In the pulse Doppler radar signal processor, when 
estimating the "average" windspeed in a given range cell, 
there may be a question as to whether "average" should be the 
statistical mean or the statistical mode (most probable 
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value) . The pulse pair algorithm will estimate the 
statistical mean. In the symmetric case the mean and the mode 
are the same. However, with skewness in the spectrum this is 
not true as seen in Figure 29 which shows that these two 
values can differ very largely due to the increased skewness 
and spectrum width. Figure 30 illustrates the difference 
between the mode of the skewed spectrum and the pulse pair 
mean as a function of skew and spectrum width. As the 
spectrum width increases, this difference is more sensitive 
to spectrum skew. Therefore, the pulse pair mean estimate 
tends not to be a good mode estimator for broad spectra. The 
idea of a mode estimator is developed further in Chapter V. A 
new approach of characterizing a summary statistic of 
windspeed within a range cell is presented using a classical 
harmonic decomposition technique. This indicates potential 
for overcoming the biased mean estimation problem with a 
skewed spectrum. 
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Figure 29. Normalized Difference Value between True 
Ilean and Node of Skewed Spectrum 
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Figure 30. Normalized Difference Value between 
Pulse Pair Mean Estimate and Mode of Skewed Spectrum 





CHAPTER IV 


THE MEASUREMENT OF WINDSPEED GRADIENT 
WITH LOW PRF RADAR 

Introduction 

A satisfactory airborne radar sensor for the detection 
of windshear must possess good ranging capability, both in 
terms of maximum useable range and range resolution. 
Otherwise the look ahead range will be too short to provide 
adequate lead time warning of a windshear condition or a 
hazardous windshear will not be resolved. This suggests a low 
to medium PRF radar. On the other hand, a capability for 
unambiguously measuring large Doppler frequencies associated 
with high windspeed conditions would normally require a 
higher PRF radar . With adequate range resolution capability, 
one approach to estimating the change in windspeed throughout 
a region is to use the difference in mean Doppler from range 
cell to range cell as a first order estimate of the spatial 
wind velocity gradient. With this approach, it may not be 
necessary to estimate large absolute windspeeds in any 
particular range cell, but only the difference in average 
windspeeds between range cells, i.e., the windspeed gradient. 
The magnitude of the windspeed gradient determined from 
adjacent range cells should be substantially less than the 
absolute mean windspeed in any given range cell, even for 
high turbulence environments. 
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Thus, it appears that it may be possible to relax the 
PRF requirement for large unambiguous Doppler measurement, 
allowing a lower PRF to be used to measure the windspeed 
. This would inherently allow good ranging capability 
without compromising the value of the system for detecting 
windshear . The discussion which follows develops this premise 
by illustrating the effect of measuring Doppler difference 
frequency with a reduced PRF. Simulated I and Q data are used 
to investigate the validity of this technique. It is 
anticipated that many hazardous windshear conditions can be 
detected with a reduced PRF if this technique is indeed 
proved to be valid. Relatively simple modification of the 
signal processing associated with conventional weather radar 
operating at reduced PRF could enable implementation of the 
technique . 


Aliasing Effect 

With a given PRF, unambiguous Doppler measurement is 
possible within a ±0.5 PRF range of frequencies. In 
meteorological measurements where a much wider range of 
Doppler return frequencies is likely when high windspeed 
conditions exist, the radar return Doppler which exceeds the 
PRF/2 magnitude is aliased as a lower frequency. Thus, for 
highly turbulent conditions where the spread of Doppler 
frequencies is very high and significant aliasing can take 
place, the ratio of the peak magnitude of the return at the 
"mean" Doppler frequency to the level of the return at other 
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frequencies may be significantly degraded by aliasing. 
Obviously, reducing the PRF would further degrade this 
"signal-to-noise" ratio and may reduce the ability to locate 
the peak magnitude of the estimated spectrum of the return. 
Of course a reduced PRF directly reduces the signal-to-noise 
ratio in a given range cell simply because for a fixed 
observation time the number of signal returns is reduced. 

Of perhaps even more significance is the effect of 
aliasing on the measurement of spectral width of the Doppler 
return. If the reduced PRF Nyquist interval is inadequate to 
unambiguously represent the breadth of Doppler frequencies 
present, there can be no way of estimating this width. 
Coupled with a reduced signal-to-noise ratio, the ability to 
detect turbulent windspeed conditions within a range cell 
with reduced PRF radar appears futile. However, the utility 
of low PRF in measuring large scale windspeed gradients 
between range cells may not be so limited. 

To evaluate the effect of aliasing, two different 
scenarios will be examined. In each situation an assumed 
return with a postulated spectrum will be represented as an 
I/Q data sequence. This sequence will be examined with DFT 
processing at a rate associated with a high PRF weather radar 
processor. The sequence will then be decimated and 
reprocessed at a rate which would correspond to a reduced 
PRF. The decimation process not only reduces the number of 
data samples, but also produces the same data that would have 
been obtained had a reduced PRF radar been employed. 
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The first situation includes analysis of a pure complex 
sinusoid with no added noise. The complex sinusoid represents 
a specular signal with a one-sided bandwidth. Aliasing with a 
sinusoid signal does not affect the signal-to-noise ratio but 
only the frequency estimate of the sinusoid itself. The 
second situation examined involves an analysis of a simulated 
weather signal which has a Gaussian frequency spectrum 
typical of weather returns. This particular type signal is 
degraded by aliasing quite dramatically when considering the 
reduced PRF system. However, a useful range of mean Doppler 
difference can be preserved within the available sampling 
bandwidth . 

Analysis of a Pure Sinusoid 

Consider a situation corresponding to a noiseless 
specular pulse Doppler radar return present in two adjacent 
range cells. Fourier analysis of a record of sampled values 
will yield a spectral estimate which is corrupted by spectral 
leakage associated with the limited observation record and 
aliasing associated with the limited sampling rate [39] . In 
the radar situation the data record is limited by the 
observation time within a range cell and the sampling rate is 
just the PRF . Using the complex Discrete Fourier Transform 
(DFT) , two complex sinusoids with a difference frequency less 
than the unambiguous processing frequency range (Nyquist 
bandwidth) have been analyzed. Using a sampling frequency 
(PRF) of 250 Hz, Figure 31 illustrates the difference in the 
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Figure 31. The Frequency Difference of Aliased Complex 
Sinusoids between 500 Hz and 425 Hz 
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magnitude of a DFT processed sinusoid at 500 Hz which aliased 
to 0 Hz and a 425 Hz sinusoid which aliased to -75 Hz. Even 
though both of these signals are undersampled and the 
original frequencies are not preserved, the difference 
frequency between the two sinusoids is clearly evident . 

Theoretically, the aliasing effect on a pure sinusoid 
does not deteriorate the signal— to— noise ratio. It only 
relocates the peak and prevents unambiguous identification of 
the actual frequency. However, with spectral leakage 
associated with the DFT, there is some deterioration of 
signal-to-noise ratio in the sense that the leakage power is 
aliased within the processing bandwidth, resulting in a 
reduced peak to leakage power ratio. This leakage power is 
effectively cancelled out in the difference processing 
associated with Figure 31. As a result, any frequency 
difference which is less in magnitude than one-half the PRF 
is preserved and can be detected and estimated as the 
frequency difference between the positive peak and the 
negative peak. 

Based upon this analysis it appears that if the radar 
return from adjacent range cells is specular the Doppler 
difference is preserved and easily recognized even though an 
absolute Doppler measurement is lost through aliasing. Of 
course, these signals are not representative of weather radar 
returns in turbulent wind environments. These signals have 
been used simply to illustrate the concept. The next section 


v 
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begins a discussion of this approach using more 
representative Doppler spectra. 

Analysis of a Simulated We ather Signal 

Power Spectrum Simulation and 
Generation of I and 0 Data 

A Gaussian power spectrum can be represented by 


G 


k 


7 — e xp [ 

V 2 k G 


< f k - f m > 2 

2G 2 


] . 


Gfc is a spectral coefficient corresponding to f k ; f m and G are 
the desired mean frequency and standard deviation 
respectively. The frequency dependent signal power density 
can now be defined as S k =C-G k where C is a signal-to-noise 

power scaling constant given by 


C = 


(SNR/10) 

10 



where SNR is the signal-to-noise ratio in dB . Defining the 
receiver's white noise power per discrete frequency, N k , as 

the reciprocal of the number of spectral coefficients, 
i . e . (ZN k =l) , the power spectrum can be given by 

P k -ln(X k ).(S k+ N k ) 
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where X k is a uniformly distributed random variable over the 

interval (0,1) . For this simulation, normalization is 
performed in the frequency domain with Pjcn =p k/ p t • In order to 

transform the power spectrum into the time domain, the real 
and imaginary components of spectrum were obtained from 

V^kn 003 * 271 ^)' B k = P kn Sin < 27t y,<> 

where y k is a uniformly distributed random variable over the 

interval (0, 1) . From these, we can obtain the in-phase and 
quadrature components through 

1 (i) + jQ (i) ( A^ + jB^ ) exp ( j2jlki/m) 

k = l 


where m should be larger than 30 to avoid any aliasing in the 
time domain [40]. 

In this particular spectral model the peak of the 
spectrum represents the strongest Doppler return from within 
the range cell and is interpreted as the mean windspeed 
within that cell. The variance of the Doppler spectrum is a 
measure of the distribution of windspeeds within the cell and 
will be larger in the presence of turbulence. A less 
turbulent condition would have a spectrum which is more 
peaked and which more closely resembles the specular 
situation discussed in the previous section. 
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The investigation here is concerned with the ability to 
measure the difference between spectral peaks when these 
modelled spectra are subjected to aliasing and differenced. 
The next subsection analyzes this further. 

Typical Weather Radar Signals 

A data record consisting of 1024 sample points of a 
representative power spectrum were simulated and I,Q data 
were obtained from that power spectrum. The simulated signals 
nre assumed to be obtained with a high PRF radar without any 
serious signal aliasing problems. To analyze the effect of 
low PRF, one fourth of the simulated signal data were taken 
(which means one fourth of the original PRF) and processed 
using the DFT in the same way as the complex sinusoids in 
Figure 31. Figure 32 is a processed simulated typical weather 
radar return signal corresponding to a PRF of 2000 Hz and 
Figures 33 through 35 are aliased versions of this type 
signal as might be obtained with a low PRF of 500 Hz. A 
different mean frequency is associated with each of these 
Figures varying from 300 to 475 Hz to illustrate the effects 
of peak relocation by aliasing and reduced peak to background 
ratio associated with spectral leakage. Various spectrum 
magnitude differences are plotted in Figures 36 through 38 
with a Nyquist interval associated with the reduced PRF. In 
each case the aliased spectral peaks are evident in the 
difference with the background level between peaks 
effectively cancelled. The original difference between 
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Figure 34 . The Aliased Weather Signal with the Mean 
Frequency of +400 Hz (PRF=500 Hz) 






35. The Aliased Weather Signal with the Mean 
cy of +475 Hz (PRF=500 Hz) 
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spectral mean values is preserved and can be observed as the 
frequency difference between the maximum peak magnitude and 
the minimum peak magnitude on these plots. In Figure 36 this 
^ f f erence appears to be between a maximum peak at 
approximately -25 Hz and a minimum peak at approximately -200 
Hz yielding a difference of 175 Hz which is the difference 
between the peaks of the unaliased spectra at 475 Hz and 300 
Hz. The ambiguity problem arises in Figure 37 in that the 
Nyquist bandwidth is inadequate to represent the mean 
difference frequency and this figure shows the peaks of the 
aliased spectra at 100 Hz and -100 Hz corresponding to the 
actual peaks at 100 Hz and 400 Hz. The true difference 
frequency magnitude is larger than PRF/2, thus it is wrongly 
interpreted as a 200 Hz difference. A true frequency 
difference of 200 Hz between two aliased peaks at -100 Hz and 
200 Hz is shown in Figure 38. 

These examples are intended to represent computed 
differences in the spectral estimate of one range cell and 
that of an adjacent range cell. It appears that the detection 
of mean windspeed changes between adjoining range cells 
through low PRF radar is quite possible though signal peaks 
are not as well defined as in the case of noiseless signals. 

Discussions and Problems 

Based upon this analysis it can be readily seen that 
reduced PRF introduces aliasing and the peak to background 
ratio (SNR) will decrease accordingly, thus circumventing the 
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Figure 36. The Frequency Difference of Aliased Weather 
Signals between 475 Hz and 300 Hz 
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Figure 37. The Frequency Difference of Aliased Weathe 
Signals between 400 Hz and 100 Hz 
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Figure 38. The Frequency Difference of Aliased Weather 
Signals between 400 Hz and 200 Hz 
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accurate estimation of mean frequency and spectrum width. 
Comparison of Figures 32 and 33 shows that the simulated 
weather signal of a low PRF radar is more contaminated with 
the background noise due to aliasing of high frequency 
components . 

In spite of this aliasing problem, it appears that the 
mean frequency difference between aliased spectra can be 
preserved and detected because the peak location of the 
Gaussian spectrum is changed in a predictable way and the 
broader band incoherent background power level tends to 
cancel out. It can be noted that a difference in mean 
frequency larger than PRF/2 causes ambiguities. This 
demonstration suggests that the DFT magnitude difference of 
low PRF radar signals between range cells can overcome the 
difficulties of decreased SNR by aliasing and present the 
velocity gradient rather clearly. 

In evaluating the utility of pulse Doppler weather radar 
for detection of windshear from an airborne platform the 
motivation for low PRF was stated earlier as a desire for 
large unambiguous range capability. Even though the analysis 
presented here tends to support the premise that windspeed 
gradient can be estimated within the Nyquist bandwidth 
associated with low PRF, other problems may be predominant. 
Turbulent wind conditions within a range resolution cell can 
cause a radar return which has a quite broad spectrum with a 
very small mode. Figure 39 illustrates the effect of aliasing 
a simulated weather spectrum which does not have a well 
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Figure 39, The Aliased Signal Obtained with PRF=250 Hz 
(mean freq.=300 Hz, S/N=30 dB, width=18.75 Hz) 
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recognized peak. The peak to background level ratio is so low 
that the mean of the spectrum can no longer be identified in 
Figure 39. As shown in Figure 40 the difference between two 
such spectra yields a result which appears useless for 
estimating windspeed gradient. 

Another consideration with low PRF radar is that the 
total magnitude of the return may be small in a turbulent 
situation and any reduction in return signal level because of 
a limited number of returns in a processing interval will 
cause the signal-to-noise ratio to be much too low for any 
practical use. Finally, a problem with the use of radar on an 
airborne platform in the neighborhood of urban airports is 
the clutter environment. High clutter levels may create a 
signal to clutter level ratio which precludes the use of a 
low PRF. Therefore, even though the analysis presented in 
this chapter supports the idea that a reduced Nyquist 
interval may be adequate for estimating windspeed gradients 
associated with windshear, it should not be concluded that 
simply reducing the weather radar PRF is the best means of 
improving the ranging capability of a weather radar. 





CHAPTER V 


NEW APPROACH IN THE DETECTION OF 
HAZARDOUS WINDSHEAR CONDITIONS 

Introduction 

The aviational hazard often caused by microbursts can 
frequently be identified from a Doppler radar return by an S 
curve characteristic which describes mean windspeed changes 
along the radar range radial. The mean value of the weather 
return spectrum is generally considered as representing the 
windspeed in each range cell. However, based upon the results 
from Chapter III, in the skewed spectrum case or with the 
multimodal return spectrum, the modes of spectrum may provide 
more reliable information than the statistical mean for the 
purpose of windshear detection. Therefore, the mode 
estimation technique using the modified Prony method is 
presented in this Chapter. Also this mode estimator may be 
useful in recognizing the windshear hazard without the need 
for clutter rejection filtering. 

As has been noted earlier, one of the more popular 
methods of estimating mean Doppler or mean wind speed within 
a range resolution cell is the pulse pair estimation 
technique. This is computationally much more efficient than 
DFT based methods although spectrum parameter estimates 
involving the DFT are generally considered to be more robust. 
With an airborne Doppler radar wet microburst return, where 
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signal-to-clutter ratio is large enough, the simple pulse 
pair estimator, for example, generally yields a very accurate 
estimate of mean wind velocity in each range cell. However, 
difficulties arise in a dry microburst case since a very low 
signal-to-clutter ratio may seriously bias the mean velocity 
estimates without effective and efficient clutter filtering. 
To make matters worse, the removal of clutter may not be an 
easy task though several methods have proven to be useful 
[ 4 5 ] , [ 18 ] , [ 19 ] . It has been shown [45] that efficient clutter 
suppression can be done using an auto-regressive least 
squares method, but mean estimates from clutter-only range 
cells often fluctuate randomly, because there remain only 
weak background noise signal? after filtration. This can also 
occur when the weather return spectrum falls largely within 
the clutter filter notch and is mostly removed with clutter 
rejection processing. Another problem is that radar system 
phase noise may limit the clutter rejection capability 
yielding a too low signal-to-clutter ratio in the filtered 
spectrum thus causing an inaccurate estimation of the mean 
velocity as explained in Chapter II. 

An alternate approach to identifying the presence of a 
weather return is to locate strong peak points in the 
spectrum which may well represent the velocity spectrum modes 
of wind and clutter signals in each range cell. These peak 
values may be adequate to identify the microburst S curve 
signature and detect a hazardous windshear condition. This 
new approach is particularly attractive since it does not 
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require processing to estimate the entire spectrum but only 
involves finding a few peak points of the spectrum. The 
modified Prony method [46] is investigated here to find 
strong peak points of simulated weather spectra that include 
microburst and static clutter signals. 

Modifi ed Prony Method 

The modified Prony method [46] involves approximating a 
complex data sequence by a model consisting of undamped 
complex sinusoids. It is similar to Pisarenko Harmonic 
Decomposition (PHD) method [47], but the Prony algorithm is 
generally better than PHD procedure since it needs neither 
autocorrelation lags nor a more computationally complex eigen 
equation solution. The Prony method requires only the 
solution of two sets of simultaneous linear equations and a 
polynomial rooting. It is summarized briefly in the 
following : 

1. Find the coefficients of a complex polynomial 
minimizing the squared smoothing error. 

2. Root a complex polynomial to determine frequencies. 

3. Solve for the amplitude of each frequency. 


The 2p component Prony model is represented as 


fcau-^h ^" 1 


k=l 


(5.1) 
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where h k =A k exp ( je k ) and Z k = e xp ( j 27tf k T) . The polynomial 
constructed with roots that are the Z k of Equation (5.1) has 
the form 


2p 

*< z) =n< z -v 

k-1 


xf 1 2p-k 

=2^a[k]Z 


where a[0]— 1 by definition. Due to the unit modulus property 
Z) C _1 = Z }c * f it can be shown that the conjugate property 

a [ k ] =a [ 2p ] a * [ 2p-k ] for k=0 to k=2p must exist between the 
coefficients. Therefore, the homogeneous linear difference 
equation that has Equation (5.1) as its solution is 

a [2p] x [n-p] (a [2p-k] x [n-p+k] +a[2p]a[k] x [n-p+k] ) =0 (5.2) 

k-1 


for 2p+l<n<N, where N is the number of given data points. A 
more convenient form of (5.2) can be obtained which yields 
the conjugate symmetric difference equation 



k-l 


(g 2p [k]x [n-p+k] +g* p [k] x [n-p-k] ) =0 


where g 2p [ k ] =a [2p-k ] /a [ 2p ] and g 2p *=a[k]. Since N usually 
exceeds the minimum number needed to fit a model of 2p, i.e., 
N>2p+1, the squared smoothing error given by 
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n=p+l 


(n) 


2 


is minimized based on the measured complex data samples where 


e 2p (n)=x[n] 


♦t 


(g 2p [k]x[n+k]+g* p [k]x[n-k]) . 


Setting the complex derivatives of p 2p with respect to g 2p [l] 
through g 2p [p] to zero yields 


R 2p g 2 P = 


2 P 2p 

°p , 


( 5 . 3 ) 


where 


9 2 T p =[g 2 P [p] . • • ^g 2p [i],i,g; p [i], . . .,g; p [p] ] 
and 0 p is a pxl all zero vector. R 2p can be expressed as 
* 2p =£ [x 2 * p t n J X 2 P [ n ] + JX 2 P [ n ] x» p [n ] J] 

n=2p+l 

where J is a ( 2p + l ) x ( 2p + l ) reflection matrix and H means 
complex conjugate transposition. Here x 2p [n] is defined as 
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xL [n] = [x[n],...,x [n-p+1] , x [n-p] , x [n-p-1 ] , ... ,x [n-2p] ] . 


The fast algorithm to solve the symmetric covariance normal 
Equation (5.3) was developed by Marple [48]. 

After solving for a[l] through a[2p], the roots of the 
complex polynomial, Z k , can be obtained using the polynomial 
factoring algorithm. Then, h k in Equation (5.1) is computed 
minimizing the squared error with respect to each of the h k 
parameters, i.e., 


where 


11... 

1 ~ 


“ h 



Z 1 Z 2 

Z 2p 

, h= 

l 

h 2 

, X= 

X[l] 

x[2] 

zf 1 z!?- 1 . . 

. z!" 1 

2p 


3p_ 


_x[N]_ 


Performance and Computational Complexity 
A second order Prony model was used here to find peak 
points of simulated weather spectra. This data set had been 
previously analyzed using adaptive clutter rejection 
filtering and pulse pair mean estimation [45]. Marple's 
programs [46] were slightly modified to avoid numerical ill- 
conditioning in some cases. A 512 point complex data sequence 
from each range cell was processed. Some typical DFT spectrum 
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plots are shown in Figures 41 through 44 with the Prony 
method peak estimates also indicated. As seen in Figures 
41,42 and 43, the Prony method is able to locate spectrum 
peak points. However, Figure 44 shows somewhat inaccurately 
estimated peak points because of the presence of strong 
clutter power and the closeness of weather and clutter 
spectral peaks. 

In order to check the usefulness of this new approach 
for detection of windshear, data from 40 contiguous range 
cells which include a dry microburst with clutter were 
processed and peak velocity points were plotted versus range. 
The resulting Figure 45 clearly shows :he S curve 
characteristic around the range cell 27. 

Another important consideration with any algorithm is 
computational complexity which must not prohibit real time 
processing. Some comparisons with other spectrum estimation 
methods are made in Table 1. Of course, the Prony method is 
computationally much more complicated than other classical 
spectrum estimation techniques as the model order increases, 
but as it can be seen from Table 1, the second order Prony 
model used here requires less computation than the DFT 
method. Therefore, the modified Prony method may be useful as 
a component of a windshear detection algorithm. 




* FFT Spectrum 

• Prony Method 

Figure 41. Mode Estimates Shown in Simulated 
Weather Spectrum of Range Cell 24 



’ FFT Spectrurr 
• Prony Method 

Figure 42. Mode Estimates Shown in Simulated 
Weather Spectrum of Range Cell 29 
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' FFT Spectrum 
• Prony Method 


Figure 43. Mode Estimates Shown in Simulated 
Weather Spectrum of Range Cell 30 



■ FFT Spectrum 
♦ Prony Method 


Figure 44. Mode Estimates Shown in Simulated 
Weather Spectrum of Range Cell 25 



3 


Range Cell Number 


o a stronger mode 


a weaker mode 


Figure 45. Estimated Spectrum Modes of Simulated 
Weather Data in All Range Cells 
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Table I. Comparison of Computational Complexity Where N=512 

( f(p) means that the required number of computations depends 
on the algorithm used for a polynomial rooting) 


method 

computation 

requirement 

approximate number 
of calculations for a 
reasonable model order 

FFT 

Nlog 2 N 

4700 complex adds/mults 

AR LSQ 

2NP+P 2 

10500 complex adds/mults 
for p=10 

Prony 

2NP+18P 2 +P 3 +f (P) 

2300 complex adds/mults 
for p=2 















CHAPTER VI 


CONCLUSIONS AND RECOMMENDATIONS 
FOR FUTURE WORK 

Chapter II of this dissertation has presented a method 
of evaluating the effect of pulse Doppler radar phase 
instabilities (jitter) on the estimation of return signal 
Doppler spectrum parameters using the pulse pair and DFT 
algorithms. Phase instabilities are modeled as statistically 
independent from the radar return itself. This phase jitter 
model can in fact represent a cumulative effect from various 
sources within the radar transmitter/receiver. As seen from 
the results presented, phase jitter is not expected to cause 
any bias error in the pulse pair and DFT mean estimates. 
However, it can cause an increase in the variance of the mean 
estimate and is essentially independent of the radar return 
spectrum width. On the other hand it has been confirmed there 
is significant increase in the bias of the pulse pair width 
estimate particularly with more narrow return Doppler 
spectrum widths. Some increase in the variance of the pulse 
pair width estimate with increasing phase jitter is also 
noted. It is also confirmed on the basis of these results 
that the phase jitter effect on pulse pair estimation is 
negligible if T s is small enough and most of phase noise 

spectrum power is concentrated around the carrier frequency. 
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However, if these conditions are not satisfied, any phase 
jitter can cause significant estimation errors. 

Comparing the results of DFT estimation with the pulse 
pair estimation case, it can be said that DFT estimates are 
relatively immune to phase jitter since the DFT estimate 
errors are much smaller than the comparable results using 
pulse pair estimates. This result may be anticipated because 
it is well-known that DFT estimate variances increase much 
more slowly, as the normalized spectrum width WT 3 increases, 

than the exponential increase associated with the pulse pair 
method. One example was shown using actual phase jitter data 
to demonstrate the general approach for evaluating the system 
performance with arbitrary phase jitter spectral data. 

A significant contribution of this work is the 
development of a technique whereby radar phase instability 
measurements can actually be incorporated into an analytical 
evaluation of pulse pair and DFT estimation quality by using 
appropriate numerical techniques. This precludes the 
necessity to rely on any particular notion of a closed form 
distribution for the phase jitter spectrum. Radar system 
phase stability design specifications can now be 
quantitatively determined based on the desired quality for 
spectrum parameter estimates. 

The phase jitter analysis is based upon earlier 
published results for pulse pair spectrum parameter estimates 
which, in many cases, have been derived through a series of 
simplifying assumptions. These assumptions may limit 
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application 

of the 

results to weather spectra which 

are 

narrow with 

respect 

to the Nyquist 

processing 

interval 

and 

which are 

Gaussian 

, or at least 

symmetric. 

as well 

as 

situations 

where a 

large signal-to- 

■noise ratio 

exists . 

The 


pulse pair estimator variance expressions presented here were 
derived using perturbation analysis under the assumption that 
perturbations are not excessive (for example, a second order 
expansion was used in the derivation of (2.11)) [29]. In some 
cases this is valid only if the number of samples is very 
large, which suggests that a higher order expansion may be 
required in the derivation when considering a small number of 
samples. Therefore a re-evaluation of the assumptions made in 
deriving the error expressions for the pulse pair estimators 
may be necessary to be applied in turbulent weather 
environments if this work is to be extended further. 

One of the assumptions made in evaluating the pulse pair 
estimates is that the Doppler return spectrum is symmetric. 
This may not be valid particularly in the turbulent weather 
situation. The analysis in Chapter III shows that the mean 
estimates can be seriously biased due to skewness in the 
weather spectrum as the spectrum is broadened even though the 
width bias error can be considered to be negligible. 
Degradation of estimation quality due to the bias term is 
less than 15% if WT S is not larger than 0.15 as seen in 

Figure 27, but this condition may not always be satisfied. In 
the skewed spectrum case, the suggested poly-pulse pair 
method was demonstrated as useful in reducing bias errors of 
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mean estimates. It is also shown in Chapter III that the mode 
of the skewed spectrum and the pulse pair mean can differ 
very largely as the spectrum width increases . This may be a 
problem in some applications such as windshear detection 
where frequently the mode of a return spectrum may be a more 
informative statistic than the mean value. 

It is anticipated that windshear can be detected by 
measurement of the windspeed gradient. From the results 
presented in Chapter IV, it is theoretically possible that 
windspeed gradient can be measured with a low PRF radar if 
the difference of mean windspeed between range cells is not 
larger than (PRF-X)/4. Low PRF does introduce an aliasing 

effect which prevents the measurement of high windspeed, but 
it is shown that the gradient information, i.e., the mean 
frequency difference, can be preserved and detected. The 
difficulties of decreased SNR due to aliasing also can be 
overcome by using the magnitude difference of return spectra 
between range cells. Therefore, Chapter IV demonstrates that 
the low PRF radar has a potential for the direct measurement 
of windspeed gradient. However, in a practical situation, 
this method may have some limitations. A very broad weather 
return spectrum can cause meaningless results because of 
seriously aliased weather signals. Also the total return 
signal power may be very small in dry weather conditions and 
a limited number of returns associated with low PRF radar may 
cause too low an SNR to be useful. Furthermore, in an 
airborne Doppler radar, the high clutter return level is a 



119 


serious problem and may preclude the use of low PRF radar 
since efficient clutter filtering is almost essential in 
these situations for the detection of windshear phenomena. 

A new approach explained in Chapter V shows that 
windshear detection may be possible using a pattern 
recognition type technique by finding an "S" curve 
characteristic demonstrated here using the modified Prony 
method. From the results in Figure 45, it can be said that 
the very low order Prony model may make it possible to detect 
the windshear condition without any other preliminary 
processing. However, this new approach also has the 


limitation that 

some 

valuable 

weather 

information such as 

spectrum width 

can 

not 

be 

obtained 

without additional 

processing. Also 

as 

shown 

in 

Figure 44 

, where clutter and 


weather spectrum modes are very close together, 
identification of a weather return is an inherently difficult 
problem to solve. In these situations, the Prony method 
appears of limited use. Other more computationally 
complicated methods such as the Pisarenko harmonic 
decomposition method may be necessary. 

Future work may include an evaluation of detection and 
false alarm rates by applying methods suggested here to data 
obtained in planned system flight tests. The thresholding 
method with use of the hazard index F-factor is simple, but 
it may be more susceptible to the mean velocity measurement 
error than the method of recognizing an S curve 
characteristic associated with microbursts. This pattern 
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recognition type technique may help eventually to build a 
more intelligent system for reliable windshear detection 
since typical microburst S curve signature information can be 
accumulated and stored for future reference in a correlation 
based detection processor. 



APPENDICES 
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A ppendix A 

Rederivation of DFT Estimate Variances 
with Gaussian Phase Noise 

Theoretical expressions for the variances of mean and 
width frequency estimate derived by Berger and Groginsky can 
be found in [33]. From these equations, the variance 
expressions for a narrow Gaussian spectrum without 
consideration of phase noise are given in [30], but 
rederivation should be made to analyze the phase noise 
effect. Therefore the autocorrelation function with phase 
noise derived in [49] can be rewritten as a series expansion 
form to obtain the return spectrum S(f), i.e., 

n 

■ (T) ) 

* + n5 

n ! t, o 

where R<f>(T) represents -.he itocorrelat ion function of the 

phase noise model, <() With the assumption of Gaussian 

phase noise the spectrum can be easily obtained by Fourier 
transform. Ignoring negligible higher order terms, the 
spectrum can be approximated very accurately as 


R (X) = Se 


-2a 


-271 w I 
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where NT S represents a white noise power within the Nyquist 
bandwidth. By substituting this spectrum into the equations 
given in [33], the variance expressions including Gaussian 
phase noise are derived as 


1 -40 w 

var ( f ) = {e T [ — = + 


2 2 

4<rw w 


0 *w, 


MT 


r * nr 1 223/2 

V K t] 2k (w +w x ) V K 


4 2 2 

40 w w 

+ -==^ + ■ 

l~ 2 2 3/2 

V 2k (w +w ) 




6 2 2 

80 w 1 w 2 


I ' 2 2 3/2 

3 y 2 ft (w +w ) 


I 1 2 2 3/2 

v 2 ft (w 1 +w 2 ) 


-20 N 2 2 22 a 2 4^-2 1 N . 2 

+e 2 ( — ) T [w +20 w +2 (Twh — (T wl +— r ( 7 ) } and 

S s 1 233 12 S 


1 IN 1 N 2 

var ( w) = [a +a -a ]+ -(— ) [a 4 +a 5 -a ] h — -(-) a 

4 MW T 2 MW b MT b 

s s 



124 


where 


-4o 2 3 w 3 6/2 O 2 w w i 3 CJ 4 3 e/2 O 4 w w 2 

a = e [ — = 1 3 — — — — h — = w,-^L_ — - 


/— 2 25/2' r^"i' r’ 225/2 

8V 7t V 71 (W +w x ) 2V 7t V 7t 


(w +w 2 ) 


4 4 


+ 4/2 q 6 w w 3 12/2* o 6 w i w 2 

/ “* 2 2 5/2 ^ P 2 2 S/9 J > 

7 C (w +wj V 71 


2 2 5/2 

(w x +w 2 ) 


4 -4a 2 1 4a 4 4a 2 2 2 -1/2 4c 4 2 2 -1/2 

a,= w e [ — ;=; — H ;=: — 4 — (w +w x ) H — - 7 ( w +w 2 ) 


2^71 w 2/71 w x 


■/27C 


32 O 6 2 2 -1/2 

H = (W +Wj ] 


37271 


2 2 2 2 

2 -40 2 W (fl 4G 2 W W 1 4C 4 W W 2 

a = 2w e [ 1 ^ W. 4 — ; 4 


r i 1 1 / ' 2 2 3/2 1 * 2 2 3/2 

V K V 71 y 27t (w 4-w x ) y 27C (W 4 -w 2 ) 


8O 6 


2 2 

w 


2 2 


8O 6 W 1 W 2 


/“ ~ 1 2 2 3/2 

3y 2 tc (w 4-w 3 ) 


I 1 2 2 3/2 

y 271 ( w x 4- w 2 > 




- 20 i 4 a 4 <4 

a 4 = e [ 3 w 4- 6 a z w 1 4- 6 crw 2 4- 4 a*w 3 ] • 


a = w 

5 



125 


2 -20 2 * 2 -,2 . 2 4 /- 2 

a g = 2w e [ w +2a^w l +2o 4 w 2 +-jcr , w 3 


and 



1 

2 2 

320w T 

s 


2 2 

w T 


4 


24 


where 


2 2 

W =w + 

1 2 


2 2.2 

W =W +Af 

2 C 


2 2 3 

W =W + — Af 

3 2 


2 

c 



Appendix B 

NASA Simulati on Parameters 


*** * SIMULATION PARAMETERS***** 1. 

0 . 

A/C Distance to Touchdown (km) 7. 

Aircraft Velocity (kts) 150.0 

Glideslope Angle (deg) 3. 

No. of Complete Scans 1. 

Time Between Scans 5. 

Roll Attitude (deg) 0. 

Pitch Attitude (deg) 0. 

Yaw Attitude (deg) 0. 

Az Integration Range/2 (deg) 6.0 

Az Integration Increment (deg) .3 

El Integration Range/2 (deg) 4.0 

El Integration Increment (deg) .2 

Rng Integration Increment (m) 100. 

Random Number Seed (0-1) .224 

Runway Number 26. 

Right (1) or Left (2) 1. 

0 . 

* * * *MICROBURST & CLUTTER****** 0. 

0 . 

Along Track Offset from TD (km) -2. 

Cross Track Offset from TD (km) 0. 

Rain Standard Deviation (m/s) 1. 

Clutter Standard Deviation (m/s) .5 

Clutter Calc. Flag (1=ON,0=OFF) 1. 

Discrete Calc. Flag (1=ON,0=OFF) 0. 

Reflectivity Calc. Thres . (dBz) -20. 

Minimum Reflectivity (dBz) -15. 

Attenuation Code (0,1,2) 2. 

0 . 

****RADAR PARAMETERS********** 0. 

0 . 

Initial Radar Range (km) 1. 

Number of Range Cells 40. 

Antenna Az - if no scan (deg) 0. 

Azimuth Scan Range/2 (deg) 0. 

Azimuth Scan Increment (deg) 3. 

Antenna Elevation (deg) 1. 

Transmitted Power (watts) 2000. 

Frequency (GHz) 9.3 

Pulse Width (microsecs) 1. 

Pulse Interval (microsecs) 268.6 

Receiver Noise Figure (dB) 4. 

Receiver Losses (db) 3. 

Antenna Type (1 ,2 , or 3) 1. 



Antenna Radius (m) .381 

Aperture Taper Parameter .316 

RMS Trans. Phase Jitter (deg) .2 

RMS Trans. Freq Jitter (Hz) 0. 

0 . 

****SIGNAL PROCESSING********* 0. 

0 . 

Number of Pulses 512 . 

Number of A/D bits 12. 

AGC Gain Factor . 6 

Processing Threshold (dB) 4. 

Clutter Filter Code (-2 to N) 0. 

Clu'.:er Filter Cutoff (m/s) 3. 

No. of Bins for F-factor Avr . 5. 
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